/*******************************************************************************
 * Copyright 2006, 2017 Jerry Huxtable, Martin Davis
 *
 * Licensed under the Apache License, Version 2.0 (the "License");
 * you may not use this file except in compliance with the License.
 * You may obtain a copy of the License at
 *
 *     http://www.apache.org/licenses/LICENSE-2.0
 *
 * Unless required by applicable law or agreed to in writing, software
 * distributed under the License is distributed on an "AS IS" BASIS,
 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
 * See the License for the specific language governing permissions and
 * limitations under the License.
 */

/*
 * This file was semi-automatically converted from the public-domain USGS PROJ source.
 */
package org.locationtech.proj4j.proj;

import java.util.Objects;

import org.locationtech.proj4j.ProjCoordinate;
import org.locationtech.proj4j.ProjectionException;
import org.locationtech.proj4j.datum.Ellipsoid;
import org.locationtech.proj4j.util.ProjectionMath;

/**
* Oblique Mercator Projection algorithm is taken from the USGS PROJ package.
*/
public class ObliqueMercatorProjection extends CylindricalProjection {

	private final static double TOL	= 1.0e-7;

	private double lamc, lam1, phi1, lam2, phi2, Gamma, al, bl, el, singam, cosgam, sinrot, cosrot, u_0;
	private boolean ellips, rot, no_uoff;

	public ObliqueMercatorProjection() {
		ellipsoid = Ellipsoid.WGS84;
		projectionLatitude = Math.toRadians(0);
		projectionLongitude = Math.toRadians(0);
		minLongitude = Math.toRadians(-60);
		maxLongitude = Math.toRadians(60);
		minLatitude = Math.toRadians(-80);
		maxLatitude = Math.toRadians(80);
		alpha = Math.toRadians(-45);//FIXME
		initialize();
	}

	/**
	* Set up a projection suitable for State Plane Coordinates.
	*/
	public ObliqueMercatorProjection(Ellipsoid ellipsoid, double lon_0, double lat_0, double alpha, double k, double x_0, double y_0) {
		setEllipsoid(ellipsoid);
		lamc = lon_0;
		projectionLatitude = lat_0;
		this.alpha = alpha;
		scaleFactor = k;
		falseEasting = x_0;
		falseNorthing = y_0;
		initialize();
	}

	public void initialize() {
		super.initialize();
		double con, com, cosphi0, d, f, h, l, sinphi0, p, j, gamma0;

		//FIXME-setup rot, alpha, longc,lon/lat1/2
		rot = true;
    lamc = lonc;

    // true if alpha provided
    int azi = Double.isNaN(alpha) ? 0 : 1;
    // true if gamma provided
    int gzi = Double.isNaN(Gamma) ? 0 : 1;
		if (azi != 0) { // alpha specified
			if (Math.abs(alpha) <= TOL ||
				Math.abs(Math.abs(projectionLatitude) - ProjectionMath.HALFPI) <= TOL ||
				Math.abs(Math.abs(alpha) - ProjectionMath.HALFPI) <= TOL)
				throw new ProjectionException("Obl 1");
		} else {
			if (Math.abs(phi1 - phi2) <= TOL ||
				(con = Math.abs(phi1)) <= TOL ||
				Math.abs(con - ProjectionMath.HALFPI) <= TOL ||
				Math.abs(Math.abs(projectionLatitude) - ProjectionMath.HALFPI) <= TOL ||
				Math.abs(Math.abs(phi2) - ProjectionMath.HALFPI) <= TOL) throw new ProjectionException("Obl 2");
		}
		com = (spherical = es == 0.) ? 1 : Math.sqrt(one_es);
		if (Math.abs(projectionLatitude) > EPS10) {
			sinphi0 = Math.sin(projectionLatitude);
			cosphi0 = Math.cos(projectionLatitude);
			if (!spherical) {
				con = 1. - es * sinphi0 * sinphi0;
				bl = cosphi0 * cosphi0;
				bl = Math.sqrt(1. + es * bl * bl / one_es);
				al = bl * scaleFactor * com / con;
				d = bl * com / (cosphi0 * Math.sqrt(con));
			} else {
				bl = 1.;
				al = scaleFactor;
				d = 1. / cosphi0;
			}
			if ((f = d * d - 1.) <= 0.)
				f = 0.;
			else {
				f = Math.sqrt(f);
				if (projectionLatitude < 0.)
					f = -f;
			}
			el = f += d;
			if (!spherical)
				el *= Math.pow(ProjectionMath.tsfn(projectionLatitude, sinphi0, e), bl);
			else
				el *= Math.tan(.5 * (ProjectionMath.HALFPI - projectionLatitude));
		} else {
			bl = 1. / com;
			al = scaleFactor;
			el = d = f = 1.;
		}
		if (azi != 0 || gzi != 0) {
			if (azi != 0) {
				gamma0 = Math.asin(Math.sin(alpha) / d);
				if(gzi == 0) {
					Gamma = alpha;
				}
			}else {
				gamma0 = Gamma;
				alpha = Math.asin(d * Math.sin(gamma0));
			}
			projectionLongitude = lamc - Math.asin((.5 * (f - 1. / f)) *
					   Math.tan(gamma0)) / bl;
		} else {
			if (!spherical) {
				h = Math.pow(ProjectionMath.tsfn(phi1, Math.sin(phi1), e), bl);
				l = Math.pow(ProjectionMath.tsfn(phi2, Math.sin(phi2), e), bl);
			} else {
				h = Math.tan(.5 * (ProjectionMath.HALFPI - phi1));
				l = Math.tan(.5 * (ProjectionMath.HALFPI - phi2));
			}
			f = el / h;
			p = (l - h) / (l + h);
			j = el * el;
			j = (j - l * h) / (j + l * h);
			if ((con = lam1 - lam2) < -Math.PI)
				lam2 -= ProjectionMath.TWOPI;
			else if (con > Math.PI)
				lam2 += ProjectionMath.TWOPI;
			projectionLongitude = ProjectionMath.normalizeLongitude(.5 * (lam1 + lam2) - Math.atan(
			   j * Math.tan(.5 * bl * (lam1 - lam2)) / p) / bl);
			gamma0 = Math.atan(2. * Math.sin(bl * ProjectionMath.normalizeLongitude(lam1 - projectionLongitude)) /
			   (f - 1. / f));
			Gamma = Math.asin(d * Math.sin(gamma0));
			alpha = Gamma;
		}
		singam = Math.sin(gamma0);
		cosgam = Math.cos(gamma0);
		sinrot = Math.sin(Gamma);
		cosrot = Math.cos(Gamma);
		
		u_0 = no_uoff ? 0. :
			Math.abs(al * Math.atan(Math.sqrt(d * d - 1.) / cosrot) / bl);
		if (projectionLatitude < 0.)
			u_0 = - u_0;
	}

    @Override public void setGamma(double gamma) {
        this.Gamma = gamma;
    }
    
    @Override public void setNoUoff(boolean no_uoff) {
    	this.no_uoff = no_uoff;
    }

	public ProjCoordinate project(double lam, double phi, ProjCoordinate xy) {
		double con, q, s, ul, us, vl, vs;

		vl = Math.sin(bl * lam);
		if (Math.abs(Math.abs(phi) - ProjectionMath.HALFPI) <= EPS10) {
			ul = phi < 0. ? -singam : singam;
			us = al * phi / bl;
		} else {
			q = el / (!spherical ? Math.pow(ProjectionMath.tsfn(phi, Math.sin(phi), e), bl)
				: Math.tan(.5 * (ProjectionMath.HALFPI - phi)));
			s = .5 * (q - 1. / q);
			ul = 2. * (s * singam - vl * cosgam) / (q + 1. / q);
			con = Math.cos(bl * lam);
			if (Math.abs(con) >= TOL) {
				us = al * Math.atan((s * cosgam + vl * singam) / con) / bl;
				if (con < 0.)
					us += Math.PI * al / bl;
			} else
				us = al * bl * lam;
		}
		if (Math.abs(Math.abs(ul) - 1.) <= EPS10) throw new ProjectionException("Obl 3");
		vs = .5 * al * Math.log((1. - ul) / (1. + ul)) / bl;
		us -= u_0;
		if (!rot) {
			xy.x = us;
			xy.y = vs;
		} else {
			xy.x = vs * cosrot + us * sinrot;
			xy.y = us * cosrot - vs * sinrot;
		}
		return xy;
	}

	public ProjCoordinate projectInverse(double x, double y, ProjCoordinate lp) {
		double q, s, ul, us, vl, vs;

		if (! rot) {
			us = x;
			vs = y;
		} else {
			vs = x * cosrot - y * sinrot;
			us = y * cosrot + x * sinrot;
		}
		us += u_0;
		q = Math.exp(- bl * vs / al);
		s = .5 * (q - 1. / q);
		vl = Math.sin(bl * us / al);
		ul = 2. * (vl * cosgam + s * singam) / (q + 1. / q);
		if (Math.abs(Math.abs(ul) - 1.) < EPS10) {
			lp.x = 0.;
			lp.y = ul < 0. ? -ProjectionMath.HALFPI : ProjectionMath.HALFPI;
		} else {
			lp.y = el / Math.sqrt((1. + ul) / (1. - ul));
			if (!spherical) {
				lp.y = ProjectionMath.phi2(Math.pow(lp.y, 1. / bl), e);
			} else
				lp.y = ProjectionMath.HALFPI - 2. * Math.atan(lp.y);
			lp.x = - Math.atan2((s * cosgam -
				vl * singam), Math.cos(bl * us / al)) / bl;
		}
		return lp;
	}

	public boolean hasInverse() {
		return true;
	}

	public String toString() {
		return "Oblique Mercator";
	}

	@Override
	public boolean equals(Object that) {
			if (this == that) {
					return true;
			}
			if (that instanceof ObliqueMercatorProjection) {
					ObliqueMercatorProjection p = (ObliqueMercatorProjection) that;
					return (
						Gamma == p.Gamma &&
						alpha == p.alpha &&
						lonc == p.lonc &&
						super.equals(that));
			}
			return false;
	}

	@Override
	public int hashCode() {
			return Objects.hash(Gamma, alpha, lonc, super.hashCode());
	}
}
